theme_bw() +
scale_x_continuous(limits = c(3, 7.25),
breaks = c(3, 4, 5, 6, 7),
labels = c("3\n(1,000)", "4\n(10,000)", "5\n(100,000)", "6\n(1,000,000)", "7\n(10,000,000)")) +
scale_colour_grey() +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
axis.line.y.left = element_blank(),
legend.position = "right",
strip.background = element_blank(),
axis.line = element_line(colour = "black"),
panel.border = element_blank(),
axis.title.y = element_blank()) +
labs(size='Observations', size = 4)  +
xlab("Population, log scale")
brazil_plot
#AUGUST 1, 2022
#REPLICATION FILE FOR "CITY SIZE AND PUBLIC SERVICE ACCESS: EVIDENCE FROM BRAZIL AND INDONESIA"
#PERSPECTIVES ON POLITICS
#ALISON POST AND NICHOLAS KUIPERS
#load packages
library(foreign)
library(dplyr)
library(stargazer)
library(lmtest)
library(multiwayvcov)
library(sandwich)
library(stringr)
library(ggplot2)
library(tidyverse)
library(scales)
library(estimatr)
library(lm.beta)
library(standardize)
library(modelsummary)
rm(list=ls())
#set working direction
setwd("~/Dropbox/Cities project/pk_replication_file")
#load the data
brazil_data <- read.csv("./_3_data/brazil_data.csv")
indonesia_data <- read.csv("./_3_data/indonesia_data.csv")
#Main manuscript figures and tables---------------------------------------------
#Figure 1
source("./_2_code/figure_1.R")
#Main manuscript figures and tables---------------------------------------------
#Figure 1
source("./_2_code/figure_1.R")
#Main manuscript figures and tables---------------------------------------------
#Figure 1
source("./_2_code/figure_1.R")
#AUGUST 1, 2022
#REPLICATION FILE FOR "CITY SIZE AND PUBLIC SERVICE ACCESS: EVIDENCE FROM BRAZIL AND INDONESIA"
#PERSPECTIVES ON POLITICS
#ALISON POST AND NICHOLAS KUIPERS
#load packages
library(foreign)
library(dplyr)
library(stargazer)
library(lmtest)
library(multiwayvcov)
library(sandwich)
library(stringr)
library(ggplot2)
library(tidyverse)
library(scales)
library(estimatr)
library(lm.beta)
library(standardize)
library(modelsummary)
rm(list=ls())
#set working direction
setwd("~/Dropbox/Cities project/pk_replication_file")
#load the data
brazil_data <- read.csv("./_3_data/brazil_data.csv")
indonesia_data <- read.csv("./_3_data/indonesia_data.csv")
#Main manuscript figures and tables---------------------------------------------
#Figure 1
source("./_2_code/figure_1.R")
#AUGUST 1, 2022
#REPLICATION FILE FOR "CITY SIZE AND PUBLIC SERVICE ACCESS: EVIDENCE FROM BRAZIL AND INDONESIA"
#PERSPECTIVES ON POLITICS
#ALISON POST AND NICHOLAS KUIPERS
#load packages
library(foreign)
library(dplyr)
library(stargazer)
library(lmtest)
library(multiwayvcov)
library(sandwich)
library(stringr)
library(ggplot2)
library(tidyverse)
library(scales)
library(estimatr)
library(lm.beta)
library(standardize)
library(modelsummary)
rm(list=ls())
#set working direction
setwd("~/Dropbox/Cities project/pk_replication_file")
#load the data
brazil_data <- read.csv("./_3_data/brazil_data.csv")
indonesia_data <- read.csv("./_3_data/indonesia_data.csv")
brazil_bins <-
brazil_data %>%
filter(urb == 1) %>%
mutate(pop_logged = log10(pop)) %>%
mutate(bins = cut(pop_logged, breaks = seq(from = min(pop_logged), to = max(pop_logged), by = (max(pop_logged) - min(pop_logged))/20), dig.lab = 6)) %>%
group_by(bins) %>%
mutate(n_obs = n()) %>%
mutate(prop_comm_health_avg =  mean(prop_comm_health, na.rm = T),
prop_doctors_avg = mean(prop_doctos, na.rm = T),
schools_prop_avg = mean(schools_prop, na.rm = T)) %>%
dplyr::select(n_obs, prop_comm_health_avg, prop_doctors_avg, schools_prop_avg) %>%
distinct() %>%
filter(!is.na(bins)) %>%
pivot_longer(cols = prop_comm_health_avg:schools_prop_avg, names_to = "names", values_to = "values") %>%
mutate(labels = case_when(names == "prop_comm_health_avg" ~ "Health Centers per 1000 Residents",
names == "prop_doctors_avg" ~ "Doctors per 1000 Residents",
names == "schools_prop_avg" ~ "Schools per 1000 Residents",
TRUE ~ NA_character_
)) %>%
mutate(bins_split = bins,
bins_split = str_replace_all(bins_split, "\\(|\\)|\\[|\\]", "")) %>%
separate(bins_split, into = c("start_bin", "end_bin"), sep = ",") %>%
mutate_at(vars(start_bin, end_bin), funs(as.numeric(.))) %>%
mutate(bin_position = (start_bin + end_bin)/2) %>%
ungroup() %>%
dplyr::select(bin_position, labels, values, n_obs)
brazil_bins <-
brazil_data %>%
filter(urb == 1) %>%
mutate(pop_logged = log10(pop)) %>%
mutate(bins = cut(pop_logged, breaks = seq(from = min(pop_logged), to = max(pop_logged), by = (max(pop_logged) - min(pop_logged))/20), dig.lab = 6)) %>%
group_by(bins) %>%
mutate(n_obs = n()) %>%
mutate(prop_comm_health_avg =  mean(prop_comm_health, na.rm = T),
prop_doctors_avg = mean(prop_doctos, na.rm = T),
schools_prop_avg = mean(schools_prop, na.rm = T)) %>%
dplyr::select(n_obs, prop_comm_health_avg, prop_doctors_avg, schools_prop_avg) %>%
distinct() %>%
filter(!is.na(bins)) %>%
pivot_longer(cols = prop_comm_health_avg:schools_prop_avg, names_to = "names", values_to = "values") %>%
mutate(labels = case_when(names == "prop_comm_health_avg" ~ "Health Centers per 1000 Residents",
names == "prop_doctors_avg" ~ "Doctors per 1000 Residents",
names == "schools_prop_avg" ~ "Schools per 1000 Residents",
TRUE ~ NA_character_
)) %>%
mutate(bins_split = bins,
bins_split = str_replace_all(bins_split, "\\(|\\)|\\[|\\]", "")) %>%
separate(bins_split, into = c("start_bin", "end_bin"), sep = ",") %>%
mutate_at(vars(start_bin, end_bin), list(~as.numeric(.))) %>%
mutate(bin_position = (start_bin + end_bin)/2) %>%
ungroup() %>%
dplyr::select(bin_position, labels, values, n_obs)
brazil_plot <-
brazil_data %>%
filter(urb == 1) %>%
dplyr::select(pop, prop_comm_health, prop_doctos, schools_prop) %>%
pivot_longer(cols = prop_comm_health:schools_prop, names_to = "names", values_to = "values") %>%
mutate(labels = case_when(names == "prop_comm_health" ~ "Health Centers per 1000 Residents",
names == "prop_doctos" ~ "Doctors per 1000 Residents",
names == "schools_prop" ~ "Schools per 1000 Residents",
TRUE ~ NA_character_
)) %>%
filter(names != "prop_doctos") %>%
ggplot(aes(x=log10(pop), y = values)) +
#geom_point(shape = 21, fill = "lightgrey", alpha = 0.1) +
facet_wrap(labels~., scales = "free") +
geom_smooth(method = "loess", color = "black") +
geom_point(data = brazil_bins %>% filter(labels != "Doctors per 1000 Residents"), aes(x=bin_position, y = values, size = n_obs), shape = 21, alpha = 0.7, fill = "lightgrey") +
stat_summary_bin(fun='mean', bins=20,
shape = 21, fill = "lightgrey",
alpha = 0,
color='black', size=1.5, geom='point') +
theme_bw() +
scale_x_continuous(limits = c(3, 7.25),
breaks = c(3, 4, 5, 6, 7),
labels = c("3\n(1,000)", "4\n(10,000)", "5\n(100,000)", "6\n(1,000,000)", "7\n(10,000,000)")) +
scale_colour_grey() +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
axis.line.y.left = element_blank(),
legend.position = "right",
strip.background = element_blank(),
axis.line = element_line(colour = "black"),
panel.border = element_blank(),
axis.title.y = element_blank()) +
labs(size='Observations', size = 4)  +
xlab("Population, log scale")
brazil_plot
#AUGUST 1, 2022
#REPLICATION FILE FOR "CITY SIZE AND PUBLIC SERVICE ACCESS: EVIDENCE FROM BRAZIL AND INDONESIA"
#PERSPECTIVES ON POLITICS
#ALISON POST AND NICHOLAS KUIPERS
#load packages
library(foreign)
library(dplyr)
library(stargazer)
library(lmtest)
library(multiwayvcov)
library(sandwich)
library(stringr)
library(ggplot2)
library(tidyverse)
library(scales)
library(estimatr)
library(lm.beta)
library(standardize)
library(modelsummary)
rm(list=ls())
#set working direction
setwd("~/Dropbox/Cities project/pk_replication_file")
#load the data
brazil_data <- read.csv("./_3_data/brazil_data.csv")
indonesia_data <- read.csv("./_3_data/indonesia_data.csv")
brazil_bins <-
brazil_data %>%
filter(urb == 1) %>%
mutate(pop_logged = log10(pop)) %>%
mutate(bins = cut(pop_logged, breaks = seq(from = min(pop_logged), to = max(pop_logged), by = (max(pop_logged) - min(pop_logged))/20), dig.lab = 6)) %>%
group_by(bins) %>%
mutate(n_obs = n()) %>%
mutate(prop_comm_health_avg =  mean(prop_comm_health, na.rm = T),
prop_doctors_avg = mean(prop_doctos, na.rm = T),
schools_prop_avg = mean(schools_prop, na.rm = T)) %>%
dplyr::select(n_obs, prop_comm_health_avg, prop_doctors_avg, schools_prop_avg) %>%
distinct() %>%
filter(!is.na(bins)) %>%
pivot_longer(cols = prop_comm_health_avg:schools_prop_avg, names_to = "names", values_to = "values") %>%
mutate(labels = case_when(names == "prop_comm_health_avg" ~ "Health Centers per 1000 Residents",
names == "prop_doctors_avg" ~ "Doctors per 1000 Residents",
names == "schools_prop_avg" ~ "Schools per 1000 Residents",
TRUE ~ NA_character_
)) %>%
mutate(bins_split = bins,
bins_split = str_replace_all(bins_split, "\\(|\\)|\\[|\\]", "")) %>%
separate(bins_split, into = c("start_bin", "end_bin"), sep = ",") %>%
mutate_at(vars(start_bin, end_bin), list(~as.numeric(.))) %>%
mutate(bin_position = (start_bin + end_bin)/2) %>%
ungroup() %>%
dplyr::select(bin_position, labels, values, n_obs)
brazil_plot <-
brazil_data %>%
filter(urb == 1) %>%
dplyr::select(pop, prop_comm_health, prop_doctos, schools_prop) %>%
pivot_longer(cols = prop_comm_health:schools_prop, names_to = "names", values_to = "values") %>%
mutate(labels = case_when(names == "prop_comm_health" ~ "Health Centers per 1000 Residents",
names == "prop_doctos" ~ "Doctors per 1000 Residents",
names == "schools_prop" ~ "Schools per 1000 Residents",
TRUE ~ NA_character_
)) %>%
filter(names != "prop_doctos") %>%
ggplot(aes(x=log10(pop), y = values)) +
#geom_point(shape = 21, fill = "lightgrey", alpha = 0.1) +
facet_wrap(labels~., scales = "free") +
geom_smooth(method = "loess", color = "black") +
geom_point(data = brazil_bins %>% filter(labels != "Doctors per 1000 Residents"), aes(x=bin_position, y = values, size = n_obs), shape = 21, alpha = 0.7, fill = "lightgrey") +
stat_summary_bin(fun='mean', bins=20,
shape = 21, fill = "lightgrey",
alpha = 0,
color='black', size=1.5, geom='point') +
theme_bw() +
scale_x_continuous(limits = c(3, 7.25),
breaks = c(3, 4, 5, 6, 7),
labels = c("3\n(1,000)", "4\n(10,000)", "5\n(100,000)", "6\n(1,000,000)", "7\n(10,000,000)")) +
scale_colour_grey() +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
axis.line.y.left = element_blank(),
legend.position = "right",
strip.background = element_blank(),
axis.line = element_line(colour = "black"),
panel.border = element_blank(),
axis.title.y = element_blank()) +
labs(size='Observations', size = 4)  +
xlab("Population, log scale")
brazil_plot
#Figure 2
source("./_2_code/figure_2.R")
#AUGUST 1, 2022
#REPLICATION FILE FOR "CITY SIZE AND PUBLIC SERVICE ACCESS: EVIDENCE FROM BRAZIL AND INDONESIA"
#PERSPECTIVES ON POLITICS
#ALISON POST AND NICHOLAS KUIPERS
#load packages
library(foreign)
library(dplyr)
library(stargazer)
library(lmtest)
library(multiwayvcov)
library(sandwich)
library(stringr)
library(ggplot2)
library(tidyverse)
library(scales)
library(estimatr)
library(lm.beta)
library(standardize)
library(modelsummary)
rm(list=ls())
#set working direction
setwd("~/Dropbox/Cities project/pk_replication_file")
#load the data
brazil_data <- read.csv("./_3_data/brazil_data.csv")
indonesia_data <- read.csv("./_3_data/indonesia_data.csv")
#Main manuscript figures and tables---------------------------------------------
#Figure 1
source("./_2_code/figure_1.R")
#Figure 2
source("./_2_code/figure_2.R")
#Figure 3
source("./_2_code/figure_3.R")
#Figure 7
source("./_2_code/figure_7.R")
school_plot <-
read.csv("./_3_data/indonesia_school_data.csv") %>%
mutate_at(vars(num_high_school, num_priv_high, num_mid_school,num_priv_mid, num_prim_school,num_priv_prim), funs( . = ./(pop/1000))) %>%
select(pop, num_high_school_., num_priv_high_., num_mid_school_., num_priv_mid_.,  num_prim_school_.,num_priv_prim_.) %>%
pivot_longer(cols = num_high_school_.:num_prim_school_., names_to = "type", values_to = "value") %>%
mutate(type_label = case_when(str_detect(type, "priv") ~ "Private",
TRUE ~ "Public"),
level = case_when(str_detect(type, "high") ~ "High School",
str_detect(type, "mid") ~ "Middle School",
str_detect(type, "prim") ~ "Primary School")) %>%
filter(type_label == "Public") %>%
ggplot(aes(x=log10(pop), y = value)) +
geom_point(shape = 21, alpha = 0.3, fill = "lightgrey") +
geom_smooth(color = "black") +
facet_wrap(level ~ ., ncol = 3, scales = "free_y") +
stat_summary_bin(fun='mean', bins=20, fill = "lightgrey",
shape = 21, size=2.5, geom='point') +
theme_bw() +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
axis.line.y.left = element_blank(),
legend.position = "bottom",
strip.background = element_blank(),
legend.title = element_blank(),
axis.line = element_line(colour = "black"),
panel.border = element_blank()) +
ylab("Schools per 1000 Residents") +
xlab("Population (logged)")
ggsave("./_4_outputs/figure_oa_i.pdf", plot = school_plot, width = 10, height = 4)
ideb_plot <-
brazil_data %>%
filter(urb == 1) %>%
#filter(ideb_score != 0) %>%
select(pop_2010, ideb_2005) %>%
ggplot(aes(x=pop_2010, y = ideb_2005)) +
geom_smooth(method = "loess", color = "black") +
stat_summary_bin(fun='mean', bins=20,
shape = 21, fill = "lightgrey",
color='black', size=1.5, geom='point') +
theme_bw() +
scale_colour_grey() +
scale_x_continuous(trans = 'log10',
breaks = trans_breaks('log10', function(x) 10^x),
labels = trans_format('log10', math_format(10^.x))) +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
axis.line.y.left = element_blank(),
legend.position = "bottom",
strip.background = element_blank(),
legend.title = element_blank(),
axis.line = element_line(colour = "black"),
panel.border = element_blank(),
axis.title.y = element_blank()) +
xlab("Log (Population)") +
ggtitle("IDEB Score (1-10, 2011)") +
ylim(0, 10)
ggsave("./_4_outputs/figure_oa_vi.pdf", plot = ideb_plot, width = 5, height = 2.5)
indonesia_scores <-
read.csv("./_3_data/indonesia_test_scores.csv") %>%
select(-X) %>%
left_join(., indonesia_data %>% select(town_code, pop)) %>%
mutate(country = "Indonesia") %>%
ggplot(aes(x=log10(pop), y = un_score_sma)) +
geom_point(color="darkgrey", shape = 21, alpha = 0.5) +
geom_smooth(method = "lm", color = "black") +
stat_summary_bin(fun='mean', bins=20,
shape = 21, fill = "lightgrey",
color='black', size=2.5, geom='point') +
theme_bw() +
scale_colour_grey() +
scale_x_continuous(limits = c(4.5, 6.5)) +
theme(panel.grid.minor = element_blank(),
panel.grid.major.x = element_blank(),
axis.line.y.left = element_blank(),
legend.position = "bottom",
strip.background = element_blank(),
legend.title = element_blank(),
axis.line = element_line(colour = "black"),
panel.border = element_blank(),
axis.ticks.y = element_blank(),
axis.title.y = element_blank()) +
xlab("Population, log scale") +
ggtitle("National Test Scores (2015)")
ggsave("./_4_outputs/figure_oa_vii.pdf", plot = indonesia_scores, width = 5, height = 2.5)
#subset on urban municipalities
brazil_data_models <-
brazil_data[brazil_data$urb == 1,] %>%
mutate(prop_doctors = prop_doctos) %>%
mutate(elec_comp_index = elec_index)
doc_mod_brazil <- lm(prop_doctors ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
teacher_mod_brazil <- lm(teachers_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
comm_mod_brazil <- lm(prop_comm_health ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
schools_mod_brazil <- lm(schools_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
elec_mod_brazil <- lm(electr_share_2010 ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
water_mod_brazil <- lm(piped_water_share_2010 ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
mods_brazil <- list("Doctors" = doc_mod_brazil,
"Teachers" = teacher_mod_brazil,
"Health Centers" = comm_mod_brazil,
"Schools" = schools_mod_brazil,
"Electricity" = elec_mod_brazil,
"Water" = water_mod_brazil)
options(scipen = 999)
modelsummary(mods_brazil,
output = "./_4_outputs/table_oa_i.htm",
coef_omit = "Intercept|pop_growth|density|div_index|avg_income|salary|ill_rate|align",
fmt = "%.1f",
gof_map = c("nobs", "r.squared"),
coef_rename = c("log10(pop)" = "Population (logged)",
"elec_comp_index" = "Partisan Fractionalization",
"log10(pop):elec_comp_index" = "Population X Fractionalization"),
stars = T)
#subset on urban municipalities
brazil_data_models <-
brazil_data[brazil_data$urb == 1,] %>%
mutate(prop_doctors = prop_doctos) %>%
mutate(elec_comp_index = elec_index)
doc_mod_brazil <- lm(prop_doctors ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
teacher_mod_brazil <- lm(teachers_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
comm_mod_brazil <- lm(prop_comm_health ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
schools_mod_brazil <- lm(schools_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
elec_mod_brazil <- lm(electr_share_2010 ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
water_mod_brazil <- lm(piped_water_share_2010 ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
mods_brazil <- list("Doctors" = doc_mod_brazil,
"Teachers" = teacher_mod_brazil,
"Health Centers" = comm_mod_brazil,
"Schools" = schools_mod_brazil,
"Electricity" = elec_mod_brazil,
"Water" = water_mod_brazil)
options(scipen = 999)
modelsummary(mods_brazil,
output = "./_4_outputs/table_oa_i.html",
coef_omit = "Intercept|pop_growth|density|div_index|avg_income|salary|ill_rate|align",
fmt = "%.1f",
gof_map = c("nobs", "r.squared"),
coef_rename = c("log10(pop)" = "Population (logged)",
"elec_comp_index" = "Partisan Fractionalization",
"log10(pop):elec_comp_index" = "Population X Fractionalization"),
stars = T)
#subset on urban municipalities
brazil_data_models <-
brazil_data[brazil_data$urb == 1,] %>%
mutate(prop_doctors = prop_doctos) %>%
mutate(elec_comp_index = elec_index)
doc_mod_brazil <- lm(prop_doctors ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
teacher_mod_brazil <- lm(teachers_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
comm_mod_brazil <- lm(prop_comm_health ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
schools_mod_brazil <- lm(schools_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
elec_mod_brazil <- lm(electr_share_2010 ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
water_mod_brazil <- lm(piped_water_share_2010 ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_salary_mo_2010 + ill_rate_2010 + align, data = brazil_data_models)
mods_brazil <- list("Doctors" = doc_mod_brazil,
"Teachers" = teacher_mod_brazil,
"Health Centers" = comm_mod_brazil,
"Schools" = schools_mod_brazil,
"Electricity" = elec_mod_brazil,
"Water" = water_mod_brazil)
options(scipen = 999)
modelsummary(mods_brazil,
output = "./_4_outputs/table_oa_i.html",
title = "table_oa_i.html",
coef_omit = "Intercept|pop_growth|density|div_index|avg_income|salary|ill_rate|align",
fmt = "%.1f",
gof_map = c("nobs", "r.squared"),
coef_rename = c("log10(pop)" = "Population (logged)",
"elec_comp_index" = "Partisan Fractionalization",
"log10(pop):elec_comp_index" = "Population X Fractionalization"),
stars = T)
modelsummary(mods_brazil,
output = "./_4_outputs/table_oa_i.html",
title = "table_oa_i",
coef_omit = "Intercept|pop_growth|density|div_index|avg_income|salary|ill_rate|align",
fmt = "%.1f",
gof_map = c("nobs", "r.squared"),
coef_rename = c("log10(pop)" = "Population (logged)",
"elec_comp_index" = "Partisan Fractionalization",
"log10(pop):elec_comp_index" = "Population X Fractionalization"),
stars = T)
doc_mod_indo <- lm(prop_doctors ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_income, data = indonesia_data)
teacher_mod_indo <- lm(teachers_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_income, data = indonesia_data)
comm_mod_indo <- lm(prop_comm_health ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_income, data = indonesia_data)
schools_mod_indo <- lm(schools_prop ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_income, data = indonesia_data)
elec_mod_indo <- lm(prop_govt_electricity ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_income, data = indonesia_data)
water_mod_indo <- lm(prop_water ~ log10(pop)*elec_comp_index + pop_growth + density + div_index + avg_income, data = indonesia_data)
mods_indo <- list("Doctors" = doc_mod_indo,
"Teachers" = teacher_mod_indo,
"Health Centers" = comm_mod_indo,
"Schools" = schools_mod_indo,
"Electricity" = elec_mod_indo,
"Water" = water_mod_indo)
options(scipen = 999)
modelsummary(mods_indo,
output = "./_4_outputs/table_oa_ii.html",
coef_omit = "Intercept|pop_growth|density|div_index|avg_income",
fmt = "%.1f",
gof_map = c("nobs", "r.squared"),
coef_rename = c("log10(pop)" = "Population (logged)",
"elec_comp_index" = "Partisan Fractionalization",
"log10(pop):elec_comp_index" = "Population X Fractionalization"),
stars = T)
